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ABSTRACT 



This paper describes the theory and algorithms of distance transform for fuzzy subsets, called fuzzy 
distance transform (FDT). The notion of fuzzy distance is formulated by first defining the length of 
a path on a fuzzy subset and then finding the infimum of the lengths of all paths between two points. 
The length of a path n in a fuzzy subset of the n-dimensional continuous space 3?" is defined as the 
integral of fuzzy membership values along n. Generally, there are infinitely many paths between 
any two points in a fuzzy subset and it is shown that the shortest one may not exist. The fuzzy 
distance between two points is defined as the infimum of the lengths of all paths between them. 
It is demonstrated that, unlike in hard convex sets, the shortest path (when exists) between two 
points in a fuzzy convex subset is not necessarily a straight line segment. For any positive number 
9 < 1, the ^-support of a fuzzy subset is the set of all points in 3? n with membership values greater 
than or equal to 9. It is shown that, for any fuzzy subset, for any nonzero 9 < 1, fuzzy distance is 
a metric for the interior of its ^-support. It is also shown that, for any smooth fuzzy subset, fuzzy 
distance is a metric for the interior of its 0-support (referred to as 'support'). FDT is defined as 
a process on a fuzzy subset that assigns to a point its fuzzy distance from the complement of the 
support. The theoretical framework of FDT in continuous space is extended to digital cubic spaces 
and it is shown that for any fuzzy digital object, fuzzy distance is a metric for the support of the 
object. A dynamic programming based algorithm is presented for computing FDT of a fuzzy digital 
object. It is shown that the algorithm terminates in a finite number of steps and when it does so, 
it correctly computes FDT. Several potential applications of fuzzy distance transform in medical 
imaging are presented. Among these are the quantification of blood vessels and trabecular bone 
thickness in the regime of limited special resolution where these objects become fuzzy. 

Keywords: Distance transform, fuzzy subset, adjacency, path, fuzzy distance, metric, dynamic 
programming. 
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1. INTRODUCTION 

Effective tools for object shape analysis are important and useful in many imaging applications 
including the medical ones. One such popular and widely used tool is the distance transform (DT) 
[1-4] of an object. For a hard (binary) object, DT is a process that assigns a value at each location 
within the object that is simply the shortest distance between that location and the complement 
of the object. However, this notion of hard DT can not be applied on fuzzy objects [5-7] in 
a meaningful way. The notion of DT for fuzzy objects, called fuzzy distance transform (FDT), 
becomes more important in many imaging applications because we often deal with situations with 
data inaccuracies, graded object compositions, or limited image resolution (on the order of an 
object's structure size). In general, FDT will be useful, among others, in feature extraction [8], local 
thickness or scale computation [9,10], skeletonization [11-13], and morphological [14] and shape- 
based object analyses [15]. In particular, .FDT may be useful in — fault detection in integrated 
circuit chips or in computer motherboard circuits, analysis of the dynamics of hurricane etc. FDT 
will be useful in many medical imaging applications such as — (1) computation of local thickness 
of trabecular bone or of vessels, (2) morphology-based separation of anatomic structures having 
similar intensities, e.g., artery/vein separation. In this paper, we develop the theoretical and the 
algorithmic framework for FDT in both continuous and digital spaces, study its properties, present 
a dynamic programming based algorithm to compute FDT for fuzzy digital objects, and illustrate 
preliminary results. 

Over the last few decades, the notion of DT has popularly been applied for hard objects 
only. Most DT methods approximate the global Euclidean distance by propagating local distances 
between neighboring pixels. As observed in [4], this simple yet fundamental idea was first presented 
by Rosenfeld and Pfaltz [16] in 1966. In [1], they discussed about different digital distance functions 
including the city block and the chess board distances and described how they can be mixed to create 
octagonal DTs that are better approximations of the Euclidean distances. Okabe, Toriwaki, and 
Fukumara [17] investigated these issues in 3D; this was independently discovered and illustrated in 
[3]. In [2], Danielsson described how to reasonably compute the true Euclidean DT in digital spaces 
including the 3D case. Ragnemalm [18] presented Euclidean DT methods in arbitrary dimensions. 
Borgefors [3,19,15,4] extensively studied DTs for hard 3D objects. In [3], she presented many types 
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of DTs including the popular algorithm that computes DT by using different local step lengths for 
different types of neighbors. Also, she presented the integer-valued approximation of the optimal 
local step lengths. Variations of this method were presented in [20,21]. In [4], Borgefors studied the 
geometry and equations of 3D DT and presented a two pass raster scan algorithm for computing 
approximate Euclidean distance transform. The notion of gray- weighted distance was introduced 
by Rutovitz [22] in the 1960's and was used by Levi and Montanari [23] in 1970 to define a gray- 
weighted skeleton. These concepts were summarized by Rosenfeld and Kak [24] in an exercise (see 
pp. 213-214). Mordeson and Nair [25] described different ways of defining path lengths in a fuzzy 
graph. 

In this paper, we develop the theoretical framework of FDT in the n-dimensional continuous 
space JJ™. Let us consider a path ir over a fuzzy subset of K™. Unlike the case of hard sets, the 
membership values of the points through -which ir passes need to be considered to determine its 
length. More specifically, the length of 7r is the integral of fuzzy membership values along 7r. 
Generally, there are infinitely many paths in a fuzzy subset between any two points x, y e 5R™ and it 
is demonstrated that, often, it is not possible to find one shortest path. The fuzzy distance between 
two points x, y G 5R™ is defined as the infimum [26] of the lengths of all paths between them. It is 
illustrated that, unlike the case of hard convex sets, the shortest path (when exists) between two 
points in a fuzzy convex subset is not necessarily a straight line. For any positive number 6 < 1, 
the ^-support of a fuzzy subset is the set of all points in 9?" with membership values greater than 
or equal to 9. It is shown that for any fuzzy subset, for any nonzero 9 < 1, fuzzy distance is a 
metric for the interior of its ^-support. It is also shown that, for any smooth fuzzy subset, fuzzy 
distance is a metric for the interior of its 0-support (referred to as 'support'). FDT is defined as 
a process on a fuzzy subset that assigns to a point its fuzzy distance from the complement of the 
support. The theoretical framework of FDT in the continuous space is extended to cubic digital 
spaces and it is shown that for any fuzzy digital object, fuzzy distance is a metric for the support 
of the object. A dynamic programming based algorithm is presented for computing FDT of a fuzzy 
digital object. It is shown that the algorithm terminates in a finite number of steps and when it 
does so, it correctly produces the FDT image. 

This paper is organized as follows. Section 2 describes the theory and properties of FDT 



in both continuous and digital cubic spaces. The metric property of fuzzy distance is studied in 
the same section. A dynamic programming-based algorithm for computing FDT of digital objects 
is presented in Section 3. The termination of the algorithm in a finite number of steps as well as 
its correctness are established in the same section. The results of applications of the FDT method 
are presented in Section 4. Finally, the conclusions are drawn in Section 5. 

2. THEORY 

In this section, we formulate the theory and study the properties of fuzzy distance transform of 
fuzzy subsets defined either in the continuous space or in a digital space. First, we will describe the 
case for the continuous space which will then guide the formulation in a digital space. 

2.1. FDT in Continuous Space 

Let 9?" denote the n-dimensional continuous space. A fuzzy subset [5] <S of 9?" is defined as a set of 
pairs {(x, pis(x)) \ x G 9? n } where fx s ■ 3?" — >■ [0, 1] is the membership function of <S in 5R™. For any 
value 9 G [0, 1], 8-support of S, denoted by 9g(<S), is the hard subset {x \ x G 9?" and ns{x) > 9} of 
9i™. In other words, the ^-support of S is the set of all points in 9? n with membership values greater 
than or equal to 9. 0-support will often be referred to as support and be denoted by 6(<S). A fuzzy 
subset with a bounded support is called bounded. 

The following notions on fuzzy subsets are used in this paper; see [27,28] for details. A 
fuzzy subset S is a ring if ns{x) = fi{r), where r = \\x — Xq\\ for some x 0 G J?" and /} : 9? -» [0, 1] 
is a membership function. S is said to be convex if, for every three collinear points x, y, and z in 
9?" such that y lies in between x and z, fj, s (y) > min[vs(x), Ps(z)]- A convex ring is called a fuzzy 
disk. A fuzzy subset is called smooth if it is differentiable at every location x G 3? n . 

Let 5 be a (hard) subset of 9?". We shall use S to denote it complement and Interior (S) 
to denote its interior which is the largest open set contained in S. The distance transform (DT) of 
5 may be represented as an image {(x, D s (x)) | x G 9?"} on 9?" where D s is the DT value at x that 
is defined as follows. 

z> 5 (*) = inf{||s-v|| \ yeS}, (2.1) 
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where, 'inf gives the infimum of a set of positive numbers and || • || is the Euclidean norm. In the 
above equation, S should be a proper subset of 5ft". In digital images, we always deal with bounded 
objects so that S is always non empty. In the subsequent discussions, we will consider the fuzzy 
distance transform of bounded subsets only. 

Here, we define FDT of a fuzzy subset S in R n . Similar to ordinary DT, FDT is an image 
on 5R n ; we shall denote the FDT image by a set of pairs {(x,Cl s (x)) j x € $t n }, where ft s {x) is the 
FDT value at x which is defined in the following way. 

A path 7r in 5i" from a point x E R n to another point (not necessarily distinct) y G 5R n is a 
continuous function tt : [0, 1] -> 5R n such that tt(0) = x and tt(1) = y. The length of a path tt in S, 
denoted by U s (tt), is the value of the following integration 

n s («) =■ £ M^t))\^\dt. (2.2) 

Following the above equation, n 5 (7r) is the integral of membership values (in S) along tt. Note 
that, if we define the inverse path ir'(t) = tt(1 - t) of it, it can be shown that n 5 (7r) = n 5 (7r'). 
A question arises at this point on weighting the path length using membership values. Although, 
other membership-based weights, e.g., absolute derivative, may be useful in some applications, we 
use membership values themselves as weights as we intend to define distance between two points as 
the minimum material that has be to be traversed through to proceed from one point to the other. 
Thus when a path passes through a low density (low membership) region, its length increases slowly 
and the portion of the path in the complement of the support of S contributes no length. This 
approach is useful to measure regional object depth, object thickness distribution etc. 

In some applications, it may be useful to consider lengths of connected paths (paths entirely 
contained in 0(5)) only. This can be conceived by incorporating a little change in the definition of 
path length as follows 

%(*) = £&Mt))\^\dt, (2.3) 

where, 

fj, s (x), i£/jL S (x)>0, 

* ■ {2A) 
oo, otherwise. 




It can be shown that both n 5 and II5 leads to identical FDT image and Theorems 2.4 and 2.5 (see 
later in this section), originally proved II 5 , are true for U' s , also. In this paper, we use II^ to define 
fuzzy distance. 

The set of all paths from a point x G 3? n to another point y G 3?" is denoted by V(x,y). 
It may be noted that V(x,y) contains infinitely many paths. The shortest path from x G W 1 to 
y G di n in S is a path ir Xty G V{x,y) such that Us(n x , y ) < Tls(ir),Vir G V(x,y). It is worth to 
mention that there may or may not exist one shortest path between two points in a fuzzy subset 
and when it exists it may not be unique. Here, we describe a 2D example where no shortest path 
exits between two points. Let us consider a fuzzy subset of 3? 2 "with its support as a disk and let 
us pick any two points x and y near its center. Within the support of the fuzzy subset, there are 
only two membership values — a high membership value for the points on the straight line segment 
joining x and y, and a low membership value elsewhere. It is not difficult to see that the straight 
line segment xy is not the shortest path between x and y and for any other path n between the 
two points it is possible to find a shorter path by further straightening the path. Although, the 
existence of one shortest path between two points in a fuzzy subset if not guaranteed, the infimum 
of path lengths always exits and it is unique. Let Cs(x,y) denote a subset of positive real numbers 
defined as (s(x,y) = {n,s(7r)|7r G V(x,y)}, i.e., (s{x,y) is the set of all possible path lengths in 
S between x and y. The fuzzy distance from x G 9^" to y G 5?" in <S, denoted as us(x,y), is the 
infimum of C(z,y), i.e., 

u s {x,y) = infCs(x,y). (2.5) 
In the following, we show that fuzzy distance satisfies metric properties. 

Lemma 2.1 For any fuzzy subset S ofR n , fuzzy distance cj s satisfies the following properties for 
any x, y G 5R n , 

1. co s (x,y) = 0, ifx^y, 

2. u s (x,y) =w <s (y,x), 

3. u s (x, y) < u>s(x, z) + u s (z, y), for any z G W 1 . 
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Proof. To prove (1), we consider the path ir(t) = x for t £ [0, 1]. Following (2.2), the length n,s(7r) 
of the path is '0' as |^| = 0 for t e [0, 1]. Hence, (1) is proved. 

Here, we prove (2). For any path tt from x to y, there is one and only one path it' from 
y to x, where ir'(t) = ir(l - t). Following (2.2), it can be shown that n,s(7r) = n 5 (7r') and thus, 
Cs(x,y) = Cs{y,x). Therefore, following (2.2), u s (x,y) = mfQ(x,y) = mi( s (y,x) = u s (y,x). 
This completes the proof of (2). 

Let z be any point in R n . To prove (3), we shall show that for any small positive value 5, it 
is possible to construct a path from x to y whose length is less than cjs(x, z) + u) S (z, y) + 5. Since, 
u s {x,z) is the fuzzy distance between x and z in S, for any small positive number e, there exists 
a path from x to z whose length is less than u> s (x, z) + e; let ir x>z be a path between x and z such 
that Iisi^x,z) < u s {x, z) + 5/2. Similarly, we can show that there is a path ir Zty from z to y such 
that n 5 (7r Z)2/ ) < u s (z, y) + 6/2. We can always construct a path it : [0, 1] -+ 3R" between x and y by 
concatenating these two paths ir XiZ and 7r z ^ as follows 

f TT xz (2t), for 0 < t < 0.5, 

7T(t) = ^ 

[ TT Zi y{2t - 1), otherwise. 

Following (2.2), 

U s (n) = J\ s (n(t))\^\dt 

= iV,^W)l^i^ + /V,K,W)l^l^ 
= n s (7r X)Z ) + n s (7r Zi!/ ) 

< w 5 (x, z) + w s (z,7/) + (5. 

The above inequality is established for any small positive number S and therefore, following that 
n 5 (7r) € Cs(x,y), the infimum of ( s (x,y) is less than or equal to tJ S (x,z) + uj s (z,y). Hence 
ws{x, y) < u s (x, z) + oj s (z, y). Since, z can be any point in 9R n , the proof of (3) is complete. □ 

Lemma 2.2 For any fuzzy subset S ofW 1 , for any nonzero positive number 9 < 1, fuzzy distance 
u s satisfies the following property for any x,y e 5R n such that either x or y is in Interior (Q e (S)), 

u s (x,y) > 0, ifx^y. 



Proof. To prove the lemma, without loss of generality, let us assume that x € Interior{Q g {S)). 
Since, x £ Interior(Q e (S)) and x and y are two definite points, there exists a positive number 5 
such that the hyperball B x5 with center at x and radius of 5 is entirely contained in Interior (Q e (S)) 
and is disjoint from y. Therefore, any path ir xy from x to y must contain a sub path -n xS from a; 
to the boundary of B xS and the Euclidean length of ir xS is at least 5. Let t be the infimum of the 
membership values over B xS . It may be noted that there may not be the smallest membership value 
over B x s but there is always the infimum. Following the fact that B x $ C Interior (Q$(S)), t > 0. 
Therefore, following (2.2), U s (ir Xiy ) > U s (tt xS ) > 50 > 0. Since, n xjy is any path from x to y, the 
fuzzy distance oj s (x, y) > 50 > 0. D 

It is interesting to note that the above theorem is not true for 0 = 0. This can be demon- 
strated using a simple 2D example. Let us consider a fuzzy subset <S with it's support forming a 
disk centered at a point x. Let y be the point on the periphery of the disk at 0 radian angle with the 
horizontal axis emanating from x. The membership values are assigned over the disk as follows — 
all points in the disk within the angular interval [0, w) (open at the upper end) around x are assigned 
1; the points within the angular interval [tt, are assigned 1/2; and so on. It is not difficult to 
see that, in this example, the fuzzy distance between the two distinct points x £ Interior{Q{S)) 
and y is 0. 

In a smooth fuzzy subset <S, the membership values over a bounded closed region X inside 
Q(S) has a nonzero infimum. The smoothness of S and the compactness of X together imply that 
the set of membership values over X includes all its limit values. Using this fact and the proof 
strategy of Lemma 2.2, the proof of the following lemma is straight forward. 

Lemma 2.3 For any smooth fuzzy subset S o/5R™, fuzzy distance u s satisfies the following property 
for any x,y eW 1 such that either x or y is in Interior(6(S)), 

us{x,y) > 0, ifx^y. 

Theorem 2.4 For any fuzzy subset S ofR n , for any nonzero positive number 0 < I, fuzzy distance 
oj s is a metric for Interior(Q e {S)). 



Proof. Obvious from Lemmas 2.1 and 2.2. 



□ 



&{S) 




Figure 1. An example to demonstrate that the shortest paths (when exist) between two points in a 
convex fuzzy subset are hot necessarily a straight line segment. The support 0(<S) of the fuzzy disk 
is the union of different shaded regions. The points within the dark grey region have a sufficiently 
high membership value in S while the points within the light grey region have a sufficiently low 
membership value. Obviously, the shortest path between the two points x and y is contained within 
. the light grey region and is not a straight line segment. 

Theorem 2.5 For any smooth fuzzy subsets ofU n , fuzzy distance lo s is a metric for Interior (9(5)). 

Proof. Obvious from Lemmas 2.1 and 2.3. □ 

It may be interesting to observe that the shortest paths (when exist) in a fuzzy subset S 
between two points x,y e 3?" are not necessarily a straight line segment even when <S is convex. 
For an example, consider the fuzzy disk S as illustrated in Figure 1. In the figure, the support 
9(5) of the fuzzy disk is the shaded region. In this example, within the support 6(5), there are 
two possible membership values — points with the high membership value are shown as dark grey 
and those with the low membership value (within the outer annular region) are shown as light grey. 
Consider the two points x and y as illustrated in the figure. Assuming the high membership value 
sufficiently close to one and the low one close to zero, the shortest path between x, y should be 
contained within the light grey region and therefore, is not a straight line segment. 

Maintaining the same spirit as in the distance transform of hard sets, the FDT value Q s (x) 
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at a point x 6 R n is equal to the fuzzy distance between x and the closest (with respect to u s ) 
point in 6(<S). In other words, the value of fi 5 (x) is defined as follows 

Q s (x) = M{uj s (x,y) \ y £ 6{S)}. (2.6) 
2.2. FDT in Digital Cubic Spaces 

Here, we describe fuzzy distance transformation in digital cubic spaces guided by the formulation 
of the same in the continuous space 5R n as discussed in the previous sub-section. To make our 
description more precise, we will redefine the concepts in digital spaces which we have already 
defined in the continuous space. In order to avoid unnecessary heaping of notations, we will use 
the same notations used in Section 2.1, although, their exact meaning in a digital space may be 
different from that in the continuous space. 

A digital space V is an ordered pair (G, a) where G is the underlying digital grid and a 
is a binary relation on G that indicates the adjacency relationship between every two points in G. 
In general, a digital grid is a set of points in 5?" such that the inter point distances are bounded 
below and within a bounded distance of any point in 5R n , there is at least one point in the grid. 
However, most imaging systems acquire images in cubic grids and these grids are simple to describe 
as well as to understand. Therefore, we confine ourselves to cubic grids only. An n-dimensional 
cubic grid may be constructed by dividing 5R n into hyper cuboids with n orthogonal families each of 
equally spaced parallel hyperplanes. The set of the centers of these hyper cuboids generates a cubic 
grid and it is not difficult to see that, under a proper coordinate system, these points represent the 
points in Z n where Z is the set of all integers. We use Z n to represent an n-dimensional cubic grid. 
The notion of adjacency is useful to define a path in a digital space and the boundary of a digital 
object. Although, researchers have used fuzzy adjacencies [29,30], the interpretation of fuzziness of 
adjacencies in the context of a path is not clear. Therefore, we confine ourselves to hard adjacencies 
only. In other words, a : Z n x Z n -> {0, 1}. Two points p,q e Z n are called adjacent if and only 
if a(p, q) = 1. In the rest of this paper, a can be standard 4- or 8-adjacency in 2D, 6-, 18-, or 
26-adjacency in 3D, and their higher dimensional analogs. Two adjacent points are often referred 
to as neighboring points to each other. 

A digital object O is a fuzzy subset {{p,Ho{p)) \ p S Z n }, where \i 0 : Z n [0, 1] specifies 
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the membership value at each point in the object. It should be noted that, in general, an imaging 
system [31] acquires images containing information of a target object, often, along with other co- 
objects. However, following the spirit of the problem addressed here, we start with the definition of 
digital objects instead of digital images [24]. In other words, the proposed framework assumes that 
the target object is already (fuzzily) segmented from an acquired digital image using an appropriate 
segmentation method [24,32-34]. This process of extracting an object from an image, commonly 
referred to as image segmentation, has extensively been studied for decades. Since the subject 
itself is an open area of research, we prefer not to enter it here, but continue describing our theory 
assuming the above definition of a digital object. The support 6(0) of a digital object O is the 
set of all points in Z n each having a nonzero object membership value, i.e., 0(0) = {p|p G 
Z n and fJ,o(p) # 0}. A path ir in a set S of points from a point p G S to another (not necessarily 
distinct) point q G S is a sequence of points (p = p\,P2, ■ ■ ■ ,Pm = <?} such that pi G S for all 
1 < i < m and pj is adjacent to p j+l for all 1 < j < m. The length of the path is m. Although, only 
hard adjacency relations are considered to define a path, we believe that the Euclidean distance 
between two adjacent points should be utilized in a meaningful way in defining its length. A set of 
points S will be called path- connected if and only if, for every two points p,q G S, there is a path 
in S from p to q. V(jp,q) will denote the set of all paths in Z n from a point p G Z n to another 
point q G Z n . For the purpose of defining the length of a path, we use the notion of a link and its 
length. A link is a path consisting of two points. The length of a link (p, q) in O may be defined 
in different ways, e.g., (1) max{/i 0 (p),/i 0 (?)} x \\p - q\\, (2) § + /xo(fl)) x \\p ~ q\\ etc. It 
may be noted that in both examples, the length of a link has two components - one coming from 
the membership values at p and at q and one other from the distance between the two points. In 
this paper, we will follow the second choice for the length of a link. Theoretical requirements for a 
valid length function for links are — (1) the length of a link (p,p) is '0', (2) the length of the links 
(p, q) and (q, p)\p,qe Z n are equal, and (3) the length of a link (p, q), where p G 6(0) and p # g, 
is greater than zero. The length II 0 (7r) of a path ir = (p = Pi,p2, • ■ ■ ,p m = q) is the sum of the 
lengths of all links on the path, i.e., 

n c W = E%W + wM x \\Pi-Pi+M ( 2 - 7 ) 
»=i 1 

n p , q G V(p,q) is one shortest path from p G Z n to q G Z n in O if Ii 0 {^ P , q ) < n 0 (7r),V7r G V(p,q). 
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Unlike the case of the continuous case, one shortest path always exists between two points in a 
bounded digital fuzzy object. 

Proposition 2.6 For any digital space V = (Z n ,a), for any fuzzy object O onV with a bounded 
support, for any two points p,q G Z n , there exists one shortest path from p to q. 

Proof. Let us consider a path ir = (pi,p 2 , ■ ■ ■ ,Pi,Ph,Pi 2 ,- ■ ■ ,Pi t ,Pi,Pi+i> ■ ' ' >Pm) with a repetition 
at pi forming a loop (a path with the common starting and finishing point). Obviously, the length 
of the path (pi,P2, 1 ■ ■ ,Pi,Pi+i, 4 " • ,Pm), obtained by removing the repetition, is less than or equal 
to that of 7T. Therefore, to prove this proposition, we show that among the paths without any 
repetition, there is one with the smallest length. Since the support of O is bounded, there are only 
finitely many links with nonzero values and therefore, there is finitely many length values for the 
paths without any repetition. Hence, there is the minimum length for the paths between p and q 
and a path with the minimum length is one shortest path. □ 

Although, the existence of one shortest path is guaranteed in fuzzy digital objects, as in 
the continuous space, there may be multiple shortest paths between two points in a fuzzy digital 
object between two points. In the rest of this paper, we consider fuzzy digital objects with bounded 
supports. The fuzzy distance, from p € Z n to q € Z n in O, denoted as u 0 (p, q), is the length of any 
shortest path in O from p to q. Therefore, 

uo{p,q)= min il 0 (7r). (2.8) 

In the following, we show that, similar to the case of fuzzy subsets in the continuous space, 
fuzzy distance in a digital object defined over a digital cubic space is a metric. 

Lemma 2.7 For any digital cubic space V = (Z n ,a), for any digital object O, fuzzy distance u 0 
satisfies the following properties for any p,q e Z n , 

1. uj 0 (p,q) = o, ifp = q, 

2. u 0 (p,q) = u 0 (q,p), 
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3. ojo(p, <?) < w 0 (p, r) + u 0 (r, q), for any r G Z n . 



Proof. Strategy of the proof for this lemma is similar to that of Lemma 2.1. To prove (1), we 
consider the path (p,p) from p to p itself. Following (2.7), the length U 0 ((p,p)) — p-oip) x = 
0. Hence, (1) is proved. 

We prove (2) as follows. For any path ir = (p = Pi,p2, ■ ■ • )P m = ?) from p to 9, we can 
construct one and only one path ir' from q to p as follows: V = (q — p m ,p m -i, • • • iPi = p)- It is 
not difficult to see that, following (2.7), II 0 (7r) = rio(7r'). Therefore, following (2.8), u>o(p,q) = 
min^ e p( P)g ) Uo(tt) = mhvePfa,?) 11© (ff') = uj 0 {q,p). This completes the proof of (2). 

Here, we prove (3) by considering any point r in Z n and assuming that 7r P)7 . = (p — 
PiiP2i • ■ • iPmp — J") is one shortest path from p to r while ir T<q = (r = 51, q 2 , ■ ■ ■ , q mq = q) is one 
shortest path from r to g. Since, p mp = <?i = r, we can concatenate these two paths 7r PiJ . and 7r r>9 to 
construct a new path ir = (p = pi,p 2 , ■ • ■ ,Pm p ,92, ■ • ■ ,9m, =9) from p to g. Following (2.7), 
m p -i 1 

n oW = E 2^°^ + ^°^+ 1 )) X llft-Pi+lll 
m,-l ^ 

+ E 2^(®) + W>(ft+i)) x II* - ft+ill 

= n 0 (7r p , r ) + no(7r r ,,) 

= lu 0 (p, r) + w G (r, 9). 

Obviously, 7r G "P(p, 9) so that w 0 (p, q) = mijVeP(p,?) n 0 (7r') < n 0 (7r) = u) 0 (p, r) +w 0 (r, g). Since, 
r can be any point in Z n , the proof of (3) is complete. □ 

Lemma 2.8 For any digital cubic space V = (Z n ,a), for any digital object O, fuzzy distance w c 
satisfies the following property for any p,q € Z n such that either p or q is in 0{0), 

uo(p,q) #0, ifp^q. 

Proof. The strategy of the proof for this lemma is similar to that of Lemma 2.2. Without loss of 
generality, let us assume that p e 6(C). Let 7r P)9 be a shortest path between p and q. Let r be the 
first point that follows p on 7r pw but is not identical to p. Since, p ^ q, this assumption is always 
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valid. Following (2.7), n 0 (7r,, i9 ) > \{p, 0 {p) + Vo{r)) x \\p - r\\. Since p ^ r, \\p - r|| > 0. Again, 
since p G 0(C), Ho{p) > 0 and obviously, no{r) > 0. Therefore, \{p, 0 (j>) + Mo( r )) x lb - r ll > °- 
Hence, u c (p; q) = ^o(^ P , q ) > 0. □ 

Theorem 2.9 For any digital cubic space V = (Z n , a), for any digital object O, fuzzy distance u) 0 
is a metric for Q{0). 

Proof. Obvious from Lemmas 2.7 and 2.8. □ 

Following the framework of FDT for fuzzy subsets in the continuous space, the FDT value 
at a point p e Z n in a fuzzy object O over a digital space, denoted by Qo(p), is equal to the fuzzy 
distance between p and the nearest point in 0(0), otherwise. In other words, the value of Qq(p) is 
defined as follows 

fto(p) = min u>o(p,q). (2.9) 
?ee(o) 

3. ALGORITHMS 

In this section, we present an algorithm for computing the FDT of digital objects. Also, we will 
prove that the algorithm terminates in a finite number of steps and when it does so it outputs the 
FDT of the digital object. Assuming a uniform neighborhood, we can think that, with respect to 
a point p e Z n , all its adjacent neighbors are ranked. A vector 6, called resolution vector, is used 
whose zth element gives the continuous distance between a point and its ith adjacent neighbor. 
Here, we use the (3" — l)-adjacency relation (i.e., 8-adjacency in 2D and 26-adjacency in 3D). 
Therefore, S is a (3™ — l)-dimensional vector. The information about the resolution vector 5 may 
be directly obtained from the application imaging system. Let N*(p) denote the set of adjacent 
neighbors of a point p € Z n excluding p itself. 

As demonstrated in [4] , a raster scan approach effectively computes regular DT for binary 
images using only two scans. A basic reason behind this is the fact that, in a binary image, the 
shortest path from a point to the background (the complement of the binary object) is always a 
straight line segment (in a digital sense). However, this is not true for fuzzy digital objects (see 
Figure 2). This makes a raster scan based approach inappropriate in computing FDT. Here, we 
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Figure 2. An example to demonstrate why a raster scan based approach fails to compute FDT 
for fuzzy digital objects in a fixed number of scans. The support 6(0) of the fuzzy object O is 
the union of different shaded regions. The dark (light) grey region resembles the set of points with 
sufficiently high (respectively, low) membership values in O. The path -k x (7r 2 ) is the shortest path 
in O from pi (respectively, p 2 ) to the complement of 9(C). Computation of the length of ix\ (7r 2 ), 
i.e., the FDT value at p x (respectively, at p 2 ) in O needs three (respectively, five) raster scans. 

illustrate this using an example as illustrated in Figure 2. Similar to the example in Figure 1, 
in Figure 2, dark grey indicates high membership in the object O and light grey indicates low 
membership. Let us consider two points p\ and p 2 as shown in the figure. Let 7Ti be a shortest 
path from pi to Q{0) and let 7r 2 be the same but from p 2 . It is not difficult to see that three raster 
scans are needed to compute the length of iri in O (i.e., the FDT value at pi) while it takes five 
scans to do so for ix 2 (i.e., the FDT value at p 2 ). Therefore, the number of raster scans needed to 
compute the FDT of a fuzzy object is dependent on the shape of the object as demonstrated in 
Figure 2. The iterative algorithm [24] for computing the distance transform of a hard digital object 
can be used to compute FDT of fuzzy objects. However, such approaches will be inefficient in this 
application as the shortest path from a point to the background may be quite complicated and for 
each iteration the method needs to scan the entire image. 

In the following, we present a dynamic programming-based algorithm to compute FDT of 
a fuzzy digital object. Before we do so, we discuss a basic problem associated with DT in a cubic 

16 



(a) (b) (c) (d) 

Figure 3. An example to illustrate a problem in DT-based thickness computation due to even and 
odd number of digital points across the medial axis, (a), (b) Two- and three-pixel thick vertical 
bars. Object points are shown by bullets, (c), (d) Digital bars after halving resolution along each 
coordinate axis. Here too, object points are shown by bullets. Boundary points are shown by light 
grey bullets. A zero FDT value is assigned at every boundary points. 

grid and present a solution to it. For the ease of illustration, we will describe the problem in 2D, 
although, it is valid in multi-dimensional spaces. Let us take two simple examples as shown in 
Figures 3(a) and (b) — one is a vertical bar with two pixel thickness and the other is the same 
but with three pixel thickness. The object points are shown by black bullets. The largest DT 
value along a horizontal line is one for the case of Figure 3(a) and it is two for Figure 3(b). This 
gives a thickness of two and four for Figures 3(a) and (b), respectively! This inaccuracy may have 
significant importance for applications with low thickness relative to image resolution. We present 
a solution to this problem as follows. Let O = (Z n ,fi 0 ) be the fuzzy object under consideration. 
We construct a fuzzy digital object C = (Z n , fi 0 ') from O as follows. 

Step Al. /i 0 '(2» + l, 2.7 + 1) = no(i,j) for all i,j E Z n . Let us call the points in O' with coordinates 
of the form (2i + 1, 2j + 1) as primary points and among those, let us call the points with 
nonzero membership as primary object points. 

Step A2. For a non primary point p € Z n , if p is adjacent to a primary non object point, assign 
fj, 0 ,{p) = 0, otherwise assign no'ip) = the mean of the membership of all adjacent primary 
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object points. 

Finally, the resolution vector 6 is updated accordingly and a zero FDT value is assigned at all points 
in 0(C). When we apply the above two steps on a binary digital object (i.e., : Z n -} {0, 1}) of 
Figures 3(a) and (b), we obtain the objects shown in Figures 3(c) and (d), respectively. It can be 
checked that DT based thickness provides correct results in both cases. In the following we present 
an algorithm to compute FDT. 

Algorithm Compute.FDT 

Input: O = (Z n ,fj,o), oe, and 5 as defined in Section 2. 

Auxiliary Data Structures: modified fuzzy object O' — (Z n ,fj, 0 '), an image (Z n ,fi) 

representing FDT of O', and a queue Q of points. 
Output: an image (Z n , Q) representing FDT of O'. 

1. compute O' = (Z n , fj. 0 >) from O - {Z n , fi a ) following Steps Al and A2; 

2. for all p e 0(0'), set fi(p) = 0; 

3. for all p e Q(0'), set fi(p) = MAX; 

/* MAX is a large value */ 

4. for all p e 6(C) such that N*(P) D Q{0') is non empty 

5. push p into Q; 

6. while Q is not empty do 

7. remove a point p from Q; 

8. find dist min = min 9eJ v. (p) [0( g ) + S rank{Piq) x \{no>{p) + Vo> (?))]; 

/* rank(p,q) gives the rank of q in the neighborhood of p */ 

9. ifdist mia < Q{p) then 

10. set fi(p) = dist min ; 

11. push all points q € N*(p) n Q(0') into Q; 

12. output the FDT image C; 

In the following two propositions, we prove that' the algorithm compute J^DT terminates 
in a finite number of steps and when it does so it produces the FDT image. 
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Proposition 3.1 For any fuzzy object O = (Z n ,fj, 0 ) over any digital cubic space (Z n ,a), the 
algorithm compute J?DT terminates in a finite number of steps. 



Proof. The algorithm is iterative, and at each iteration in the while-do loop, it removes exactly 
one point from queue Q. Also, only the points of 0(C) visit Q. Since the number of elements in 
0(C) is always finite, the algorithm compute JFDT fails to terminate in a finite number of steps 
only when some point p € 0(C) is modified infinitely many times. Following Steps 9 and 10, the 
value at p strictly decreases after every modification. Moreover, following Steps 8, 9, and 10, during 
every modification, p is assigned the value of the length of some path in C from p to a point in 
9(C)- First, we show that p is always set to IIo'(7r) for some 7r without any repetition or loop. If 
this is not true, without loss of generality, let p be the first point to break the rule. Therefore, the 
path 7r is of the form q u q 2 , ■ • ■ ,p, ■ ■ ■ ,p. But, this is possible only when p has priorly been assigned 
the length of the path it 1 = q lt q 2 , ■ • • ,p. Following (2.7), it is obvious that U 0 '{tt') < rio'(7r). But 
this contradicts Step 9. Now, since 0(C) is finite, there are only finitely many length values for 
repetition-free paths from 0(C) to p (see the arguments in the proof of Proposition 2.3). Hence, 
any point p G Z n could be modified only finitely many times. Hence compute JFDT terminates in 
a finite number of steps. □ 

To prove the correctness of the algorithm we need the following notation. Let T p Q^(l) 
denote the set of all paths of less than or equal to / points from p to some q G 0(C). 

Proposition 3.2 For any fuzzy object O = (Z n , p, 0 ) over any digital cubic space (Z n , a), when the 
algorithm computeJFDT terminates, its output equals the FDT image (Z n ,£2o') of O' where O' is 
generated from O using Steps Al and A2. 

Proof. Initially, the algorithm starts with the points of Q(0') (see Steps 2-4). Thus, Steps 8-10 
guarantee that whenever a point p is set to dist mia (in Step 10), there is always a path from some 
q E 0(C) to p with length dist min . Therefore, for any p e Z n , £2(p) > Q 0 i(p). 

It may be observed from the algorithm that it never increases the value of any p G Z n . We 
now prove that, for any p G Z n , 0. 0 >{p) > fi(p) by using the method of induction. Obviously, for 
any p G Z n , T p Q ^ (l) is nonempty implies that p G 0(C), and following (2.7), the length of any 
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(a) (b) (c) 

Figure 4. Application of FDT on a fuzzy object, (a) A 2D slice from a 3D CTA image of a 
patient's head, (b) A fuzzy object representing the bone structures in (a). The fuzzy object was 
constructed from the image of (a) by first thresholding the bone regions and then by blurring and 
subsequently adding noise to it. (c) FDT image of the fuzzy object in (b). As visually appears in 
this figure, the ridges of FDT values follow the medial axis of the fuzzy object. 

element of T p Q^(l) is always "0". Thus, in Step 2, it is guaranteed that the algorithm sets every 
point p E Z n to the length of the shortest path of T p e^7j(l) in O' . Now let us assume that the 
algorithm sets every point p G 0(C) to the length of the shortest path of T p ^^{i - 1), for some 
i > 1, in O'. We show that the algorithm sets every point p € Z n to the length of the shortest path 
of T pj Q^7y(i) in O' in Steps 7-11. Let tt p = (q = p 1 ,p 2 , ■ ■ ■ ,Pi = p) be the shortest in O' among all 
paths of T p Q^orjii) from a point q £ Q(0') to p. Following the hypothesis of the induction method, 
at a certain iteration the algorithm sets the point pi-i to = min„. eT rio-^) < 

U 0 '((pi,P2, • ■ • at Step 10, and Step 11 ensures that, in the same iteration, p is pushed into 

Q; and upon the removal of p from Q, it is assigned (see Steps 7 to 10) a value less than or equal 
to fi(p) = ftfa-i) + <5 ran fc(Pi, Pi -i) x Wo-{Pi) + A*c(ft-i)) < n 0 '(7r p ) (following (2.7); note that 
<5ranfc( Pi> pi_i) = \\Pi "Pi-ill)- Therefore, for any p € Z n , fie?'(p) > Hence, by the results in the 

previous paragraph, for any p e Z n , Q.o>{p) = &(p)- u 
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4. RESULTS AND DISCUSSIONS 



In this section, we present a few examples of the application of FDT methods presented in Sections 
2 and 3 and will describe a method of computing local thickness using FDT in Section 4.1. The first 
example is a fuzzy object constructed from a 2D slice of a 3D computed tomographic angiography 
(CTA) image of a patient's head. The original slice is shown in Figure 4(a). The inplane resolution 
of the image is 0.25x0.25 mm 2 . A fuzzy object was generated as follows. The bone regions were 
segmented from the slice image through interactive thresholding. The fuzzy object representing the 
bone region was constructed by (1) blurring (to simulate partial volume effects) the thresholded 
image over a Gaussian kernel of radius 5 pixels, and (2) subsequently adding a correlated, zero-mean 
Gaussian noise with standard deviation equal to 10% of the pixel bone fraction value. The final 
fuzzy object is shown in Figure 4(b). The motivation behind the construction of the fuzzy object in 
this way is simply to demonstrate the results on a relatively realistic image. Figure 4(c) shows the 
FDT image as derived from the fuzzy object of Figure 4(b). In Figure 4(c), the intensity values are 
proportional to its FDT value. As visually apparent in this figure, the ridges of FDT values follow 
the medial axis of the fuzzy object suggesting that there may be potential applications of FDT in 
computing skeletons of fuzzy objects. 

4.1. Computation of Thickness 

Thickness is a useful parameter in analyzing object shape and morphology. In this section, we 
briefly describe a method to compute thickness of a fuzzy object O. 

Thickness is computed along the medial surface of O. Skeletonization [11-13] is widely 
used to generate the medial surface representation of a digital object. Following the same trend, 
we compute thickness along the skeleton of an object. It is desirable to generate the skeleton of 
an object directly from its fuzzy representation. However, this demands significant research and 
developments and is out of the support of the present paper. Therefore, we use the skeleton of the 
support 0(C) of O to compute its thickness. For this purpose, we use the skeletonization method 
described in [13]. Let 5Jfc(0(O)) denote the skeleton of 9(C). At any point p € Sk(G(<D)), the 
thickness value is computed as twice the largest FDT value in the neighborhood of p. 
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(a) 



(b) 




Figure 5. Application of FDT-based thickness computation to an arterial tree, (a) A MIP 
rendering of a 3D sub volume taken from a 3D CTA image of a patient's head (after removing 
bones) showing a portion of carotid arterial tree, (b) A MIP projection of the fuzzily segmented 
arterial tree, (c) A MIP projection of the FDT image of the 3D image shown in (b). Mean and 
standard deviation of thickness values computed along the curve skeleton of the arterial tree mask 
are 2.74 mm and 1.8 mm, respectively. 



The first example of FDT-based thickness computation is illustrated in Figure 5. Figure 
5(a) shows a maximum intensity projection (MIP) rendition of a portion of an arterial tree in a 
3D CTA image of a subject's brain vasculature. The size of the image domain is 74x217x40 and 
the voxel size is 0.32x0.32x1.25 mm 3 . The image has been rendered after removing bone using a 
recently developed method [36] in our laboratory. The mask for the arterial tree was segmented 
from the rest of the tissue using scale-based fuzzy connectedness [10]. A membership value n A {p) 
at each location within the segmented artery mask was computed as 



M(P))> Xf(p)<rn Al 
otherwise, 



(4.1) 



where, / is the image intensity function, m A and a A are mean and standard deviation of intensity 



values within the arterial mask and G m 



> an unnormalized Gaussian function with m A and 



a A as its mean and standard deviation parameters. A zero membership value is assigned at e 
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(a) 



(b) 



Figure 6. Application of the FDT-based thickness computation method to high-resolution in 
vivo 3D MR images of the human distal radius showing trabecular bone structure, (a) A 2D slice 
taken from the raw 3D MR image. The central highlighted disk represents the cross section with 
the cylindrical ROI used for computing thickness of the bone trabeculae. (b) A 3D projection of 
the BVF image computed over the ROI. Mean and standard deviation of the thickness values of 
trabecular bone over the selected ROI are 102 fj,m and 42 pm, respectively, in good agreement with 
the known thickness of human bone trabecular which is on the order 100-150 /j,m [39]. 

location outside the artery mask. A MIP rendition of the membership image of the arterial tree 
is presented in Figure 5(b). The FDT image was computed for the 3D membership image of the 
arterial tree and the MIP rendition of the FDT image is illustrated in Figure 5(c). In order to 
compute thickness, a curve skeleton of the arterial tree was computed using the method in Saha 
et al. [13]. The mean and the standard deviation of thickness (vessel diameter) values along the 
skeleton were 2.74 mm and 1.8 mm, respectively. Arterial plaques can be detected by identifying 
sudden reduction of vessel diameter while tracing along the curve skeleton of an arterial tree. 

The second example of thickness computation, illustrated in Figure 6, is a high-resolution 
in-vivo 3D MR image of the human distal radius showing the trabecular bone network. The in- 
vivo images were acquired on a 1.5T GE clinical scanner. The image size was 512x256x32 and 
the voxel size was 137x137x350 /mi 3 . A slice from the raw image is shown in Figure 6(a). A 
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cylindrical region of interest (ROI) was chosen for analysis. The central highlighted disk in Figure 
6(a) is the cross section of the ROI with that slice. The image within the ROI was preprocessed 
by deshading and noise reduction using a histogram deconvolution method [37] to produce a bone 
volume fraction (BVF) map. The spatial resolution was enhanced to 68x68x88 /im 3 using subvoxel 
classification [38]. A 3D projection of the final BVF image is presented in Figure 6(b). The FDT 
image was computed from the resolution enhanced BVF image. Thickness values were computed 
using a surface skeleton [13] of the bone mask. The mean and the standard deviation of thickness 
values along the skeleton were 102 /im and 42 /zm, respectively. 

A further demonstration of the effectiveness of the proposed FDT-based thickness compu- 
tation method involving high-resolution micro-computed tomography (/i-CT) images is presented 
in Figure 7. The basic idea was to evaluate the FDT-based thickness distributions over four regions, 
each two of them taken from the same trabecular bone sample. Two /i-CT images of different meta- 
physeal samples of human distal radius were acquired on a SANCO Medical /i-CT 20 scanner at 22 
/im isotropic resolution. Each of the two /i-CT images was processed as follows. The raw gray-scale 
image was binarized to yield bone masks. A BVF map at each location within the bone mask was 
computed using an equation similar to (4.1). The BVF images of the two samples are illustrated 
in Figures 7(a) and (b) using fuzzy shell rendering [40] supported by the 3DVIEWNIX [41] system. 
Two representative slices of the raw /i-CT images, one from each sample, are shown in Figures 7(c) 
and (d). A pair of slices, separated by 330 /im (i.e., 15 slices apart), was selected from each of 
the two samples and 2D FDT-based thickness was computed for each BVF slice image separately. 
Figure 7(e) shows two distributions of the thickness values one for each of the two different slices 
chosen from the sample shown in Figure 7(a). The mean and standard deviations of thickness values 
in one slice were 119 /im and 59 um. respectively, and the corresponding values for the other slice in 
the same sample were 116 /im and 62 /im, respectively. Figure 7(f) presents thickness distributions 
for the two slices chosen from the other sample shown in Figure 7(b). The mean and standard 
deviations of thickness values in one slice is 103 /im and 57 um. respectively, versus values for the 
other slice in the same sample of 101 /im and 54 /im, respectively. As demonstrated in Figures 7(e) 
and (f), the thickness distributions for the slices from the same sample are quite similar while they 
are significantly different for the pair of slices from the other sample. These observations are con- 
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(a) (b) 




(c) (d) 

Figure 7. Preliminary validation of FDT-based thickness method using /x-CT data of human 
distal radius trabecular bone by studying thickness distributions over four regions, each two taken 
from the same trabecular bone sample, (a), (b) Fuzzy shell rendering of the 3D BVF image of 
trabecular bone networks for the two different specimens, (c), (d) Two representative slices of the 
raw Ai-CT image, one from each of the two samples shown in (a) and (b). (e) Two FDT-based 
thickness distributions, one from each of the two different slices (330 fim apart) chosen from the 
sample shown in (a), (f) Same as (e) but for the sample in (b). 
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firmed by the mean thickness values. These preliminary experiments show that the thickness value 
computed using the proposed methods is an intrinsic property of trabecular bone at a particular 
skeletal location and subject while these values may be quite different for corresponding locations 
in different subjects. 

5. CONCLUSION 

Theory and algorithms have been presented for fuzzy distance transforms representing an extension 
of the concept of distance transforms for hard objects to the more common fuzzy objects present in 
imaging. It has been shown that fuzzy distance is a metric for the interior of the support of a fuzzy 
object. A salient feature of fuzzy objects is the property that the shortest path between two points 
is no longer a straight line. It has also been shown that the raster scan algorithm commonly used 
for distance transforms of hard objects fails in the case of fuzzy objects. A dynamic programming 
based algorithm has been presented for computing FDT of fuzzy digital objects. It has been shown 
that the algorithm terminates in a finite number of steps and when it does so, it correctly computes 
FDT. Finally, several potential applications for fuzzy distance transforms in medical imaging are 
presented, including the quantification of the thickness of vascular structures and of trabecular bone. 
The method is likely to be of significant value for structural object analysis once the technique has 
been fully validated. 
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